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ABSTRACT 

Afterglows of gamma-ray bursts are observed to produce light curves with the flux following power law 
evolution in time. However, recent observations reveal bright flares at times on the order of minutes to days. 
One proposed explanation for these flares is the interaction of a relativistic blast wave with a circumburst 
density transition. In this paper, we model this type of interaction computationally in one and two dimensions, 
using a relativistic hydrodynamics code with adaptive mesh refinement called RAM, and analytically in one 
dimension. We simulate a blast wave traveling in a stellar wind environment that encounters a sudden change 
in density, followed by a homogeneous medium, and compute the observed radiation using a synchrotron 
model. We show that flares are not observable for an encounter with a sudden density increase, such as a 
wind termination shock, nor for an encounter with a sudden density decrease. Furthermore, by extending our 
analysis to two dimensions, we are able to resolve the spreading, collimation, and edge effects of the blast wave 
as it encounters the change in circumburst medium. In all cases considered in this paper, we find that a flare 
will not be observed for any of the density changes studied. 



1. INTRODUCTION 

Since the discovery of gamma-ray bursts (GRBs), there 
has been increasing interest in their properties and behav- 
ioral characteristics. Currently, GRBs are thought to be the 
result of either a massive star collapsing ( |Woosley|1993[|Mac~ 



Fadyen & Woosley|1999) or compact binary systems merging 



( [Eichler et al.|1989| ), launching a collimated relativistic blast 
wave into the circumburst medium. The blast wave sweeps 
up, shocks, and accelerates the circumburst electrons as the 
circumburst medium slows the blast wave itself. The shock- 
accelerated electrons produce synchrotron radiation as they 
interact wi th small scale magnetic fields behind the shock 
front (e.g., Paczynski & Rhoads'1993', 'Katz'1994*, 'Sari et al. 



1996; Wijers et al. 1997; Wijers & Galama 1999; Panaitescu 
& Ku mar 2002). This radiation creates an afterglow signal 
that can be observed for days at X-ray and optical frequencies 
and for even longer time scales at radio frequencies. 

Much research has been done to understand the blast wave 
that causes the observed afterglow signal, and the evolution of 
the ultra-relativistic early time and non-relativistic late time 
phases of the blast wave are known: they can be described 
by the self -similar Blandford-McKee (BM: Blandford & Mc-] 

?||1976|) and the Sedov-Von Neumann- Taylor (ST: S edov 
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emitted light curve. The general model of the afterglow light 
curve is a smooth curve with the slope being a function of the 
density of the surrounding medium as well as the power-law 
slope of the distribution of the accelerated electron population 
at the shock front. Recent observations, though, have shown 
this model is not always sufficient. Flares in X-ra y afterglows 
were first detected with BeppoSAX ( [Frontera eFa l. 2000 ) and 
once Swift was launched in 2004, it became clear th at after- 
glow flares ( [Burrows et al.|2005t|Nousek et al. 2006; 0'Brien| 
et al.'2006') as well as optical variability in the early stages of 
the burst ( Stanek et al. 2007) were a common occurrence. To 
explain the causes of these flares, researchers have began to 
study the interaction of a blast wa ve with complex structures 
such as wind termination shocks JDai & Lu||2002t |Ramirez-| 



Pe'er & Wij ers"2006>'Eldridg e et al |20Q6 
clum ps (Ramirez-Ruiz et al.||2005| ), 



K^ 

T9? 



Neumann||T96T| |Taylor||1950] ) 
L analytical solutions describe th 



Ruiz et al.''2005l 
Nakar & Granot] 

magnetic shocks (Yost et al. 2003), collisions with nearby 
star environ ments (Mimica & Giannios 2011), and massive 
shells (Mes ler et aL|2012| ). Others have speculated the flares 
could be caused by slower ejecta catching up with the blast 
wave, r eenergizing it at later times (e.g., [R ees & Meszaros 
il998; Kumar & Piran'2000^,'Sari & Meszaros|2000MMaxEam 



_ solutions respec- 
tively. Such analytical solutions describe the radial outflow of 
a collimated blast wave at early times, as well as the spherical 
outflow of the blast wave at late times. Yet there are currently 
no exact analytical solutions for the intermediate pha ses of 
the evolution, altho ugh there are approximations (e.g., |Pe'er| 
2Q12t [Huang et"^.||199 9 for the spherical case), describing 



the dynamics of the deceleration and spreading of the blast 
wave. The spreading of the blast wave as it de collimates has 
been treated analytically (Rhoads 19991 [Granot & Piran 2012) 
and recently modeled computatio nally (Zhang & MacFadyen 
20091 [Wygoda et al.[|2QTT| [vanEerten & MacFadyen„2012, 
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Van Eerten et al.|2012|[De Colle et al.|2012[ ). 

Along with studying the evolution of the blast wave itself, 
much research has been devoted to studying the shape of the 



& Zhan g 2009; Vlasis et al. 2011), from late en gine activity 
|Pema & MacFadyen 2010; Margutti et al. 2011 ), or from de- 
layed magnetic dissipation (Gi annios,2006 ). There are many 
hypotheses of the cause of the rebrightening in the light curve, 
and in this paper, we discuss the hypothesis of a transition in 
circumburst medium from a stellar wind to a homogeneous 
medium (e.g., an interstellar medium (ISM)), including wind 
termination shocks as well as sudden lower density regions. 
These types of density changes have been numerically studied 
in much de tail and are fr om the formation of a wind reverse 
shock (e.g., [van Marie et al. 2006; Eldridge et al. 2006). 

Pe'er and Wijers (PW: Pe'er & Wijers 2006) theorized that 
the flares were the result of the blast wave radiation interfer- 
ing with the reverse shock of the stellar wind, causing a tran- 
sition peri^djnjheli^htcurve. However, Nakar and Granot 
(NG: [Nakar & Granot|2007[ ) were unable to reproduce PW's 



results using one dimensional simulations and argued that the 
flares are not caused by a wind termination shock nor a den- 
sity jump in a uniform external medium. NG concluded that 
the origin of the flares in afterglow light curves has still yet 
to be discovered. Van Eerte n et al.| ( |2QQ9| ) reconciled this dis- 
crepancy by explaining NG was correct in stating a wind ter- 
mination shock does not cause rebrightening, although the re- 
sults of PW correctl y follow from their model assumptions. 
[Mesler et al.| ( |2012| ) revisited this idea of a wind termina- 
tion shock but instead of one shock, added higher density 
shells that can realistically be expected to occur, and claimed 
that flares were observed as the blast wave encountered these 
higher density shells. 

In this paper, we address a number of circumburst medium 
interaction scenarios not previously explored. We also ex- 
tend the analysis from one dimension to two dimensions to 
understand the effects of spreading of a collimated flow at 
the shock front which was not previously considered. Specifi- 
cally, we address the following questions in our analysis using 
high resolution hydrodynamic simulations in one and two di- 
mensions: 

• When a collimated blast wave traveling in a stellar wind 
environment encounters a wind termination shock, how 
does the size of the density jump affect the dynamics 
and resulting light curve? 

• What happens when there is a density drop instead of 
a jump? Does the blast wave speed up, and in turn, 
recollimate? Will this cause a rebrightening or flare in 
the light curve? 

• When the blast wave encounters an extreme density in- 
crease, does it immediately spread outwards from this 
high energy collision and cause flares in the light curve? 
Does this sideways spreading depend on the size of the 
density jump? 

To answer these questions, we use a numerical relativistic 
hydrodynamics (RHD) code with adaptive m esh refinement 
called RAM (Zhan g & MacFad yen 2006, 2009) for our numer- 
ical simulations. We use the radiation calculation methodol- 



ogy explained in | Van Eerten et al.|p010a|b| ) to calculate syn- 
chrotron afterglow light curves and spectra. In Section [2] we 
discuss the dynamics of the simulations which is followed by 
Sectionj^discussing the resulting light curves. We summarize 
and conclude our findings in Section H] Details of the resolu- 
tion of our simulations are explained m Appendix |A| and we 
discuss the specialized case of a spike in circumburst medium 
in Appendix [B] 

2. DYNAMICS OF BLAST WAVE ENCOUNTERS 

A key aspect in understanding the role of the circumburst 
medium in the flux emitted from a blast wave is the dynamics 
of the blast wave in that medium. In this section, we discuss 
those dynamics. 

2.1. Initial Conditions 

The initial conditions for an adiabatic blast wave ("jet") 
formed by a GRB are described by the BM solution in spher- 
ical coordinates. For this paper, we use a conic section trun- 
cated at a certain opening angle for the initial setup, instead 
of the full spherical solution for, at early times, spreading has 



not yet occurred. Using the initial conditions for a jet flow- 
ing in the radial direction with an opening angle, ^o. and a jet 
energy, Ej, we obtain the isotropic equivalent energy, Etso- 



-jet 
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(1) 



The initial radius, time, and dynamics of the interaction of 
the jet with the circumburst medium are calculated using the 
BM self similar solution. For this paper, we represent the 
circumburst medium interior to the density change using a 
power-law density profile which is described as: 
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where pref is the circumburst density at a reference radius, 
Rref, and i^ = 2 or to represent the stellar wind or ISM envi- 
ronment respectively. An ISM is a region of constant density, 
and throughout this paper we use "ISM" and "homogeneous 
medium" interchangeably. The isotropic energy of the system 
is given by: 



^npextl^t^c^ 
11 -4k 



(3) 



where 7 is the fluid Lorentz factor just behind the shock, and 
r = \/27 is the Lorentz factor of the shock itself. Substituting 
Equation |2] into Equation |3] and using the approximation that 
at ultra-relativistic speeds, Rq ^ cTq, the initial radius of the 
shock is: 
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(4) 



Using these expressions, we can solve for the initial time of 
the simulation. 
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(5) 



As the jet propagates into the circumburst medium de- 
scribed by Equation [2j the Lorentz factor and radius of the 
shock over time are well known-the radius will follow a sim- 
ilar solution to Equation [5] The subscript of "FS" is used to 
denote the front of the shock, or the forward shock. 
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The Lorentz factor evolves as 
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The initial conditions at To of our numerical simulations use 
the full BM profile, but for simplicity in our analytical so- 
lution, we assume the blast wave is a homogeneous shell of 
constant density within the shell, and the back of the shell 
being at a radius of Rhack- 



Rback = RfS 1 



1 



2(3-k)Tj,s 



(8) 



Equation [s] is derived to the leading order of 1/7^ using the 
approximation that the total swept up mass is contained within 
the shell. 
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Fig. 1. — Shocks formed from a blast wave in a stellar wind encountering 
a ISM. The fluid in region 1 is the ISM into which the blast wave is trav- 
eling. Region 2 contains all of the mass swept up after the encounter with 
the density change. The fluid in regions 3 and 4 is the swept up mass prior 
to the encounter, meaning that after the encounter, the amount of mass con- 
tained in regions 3 and 4 is constant. The fluid in region 3 is that which has 
been shocked by the reverse shock, and the fluid in region 4 is the unshocked 
material that is unaware of the encounter. 



To simulate the encounter with a change in external 
medium, we model the external density profile as a piecewise 
function: 

Pwind T ^ Ifenc 
PiSM 7 ^ lenc 



Pext ■ 



(9) 



where ^yenc is the fluid Lorentz factor at which the encounter 
with the new environment occurs. The dynamics before the 
encounter are described by the equations above, but the dy- 
namics during the encounter and after the encounter are much 
different. The simulation does not necessarily return to a BM 
self similar solution for a relativistic blast wave in a homoge- 
neous medium, for the BM solution is a slow attractor and a 
perturbation in the simulation, such as a density change, can 
cause large deviations from the BM solution ( Gruzinov 2000 ). 
Also, after the encounter, the blast wave may no longer be rel- 
ativistic which prevents the blast wave from evolving towards 
the BM solution in a homogeneous medium, and a solution 
similar to ST is expected instead. 

2.2. Analytical Solution 

To analytically model the dynamics of a blast wave in a stel- 
lar wind environment encountering a density change followed 
by a homogeneous medium, we extend the derivations done 
in PW and NG. With our extensions, the analytical solution is 
applicable for a wide range of density jumps as well as density 
drops. 

For simplicity, we model the shocks of the blast wave as 
simple homogeneous shells. During the encounter, the shock 
breaks up in to three different regions-regions 2, 3, and 4 of 
Figure [T] To get the full analytical solution of the simula- 
tion, the pressures, densities, Lorentz factors, and radii of the 
three regions (regions 2, 3, and 4, which are referenced by 
subscripts from here forward in this paper) are needed. 

Across the contact discontinuity, the pressure and fluid 
Lorentz factors are equal, meaning P2=P3, and 72 = 73 = Tcd- 
Before the reverse shock passes the "back" of the shock, the 
fluid in region 4 has no knowledge of the encounter, meaning 
the density, pressure, and Lorentz factor are calculated as if 
no encounter has occurred. The rest of the fluid has experi- 
enced the encounter, and to calculate the values analytically, a 
convenient factor to define is that of NG Equation 3,ip, which 
can be derived from conservation of mass. 
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Using this factor, the Lorentz factor of the fluid behind the 
forward shock is described by a modified version of Equation 
6 from NG: 



72 = lei 




-1/2 



(12) 



where Equation [12] differs from NG Equation 6 in that 72 ex 
jenc as opposed to 72 oc 74 of NG Equation 6. 

To calculate the radius of the forward shock to the leading 
order of 1/7^, we integrate the velocity of the forward shock: 
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The radius of the contact discontinuity is calculated in the 
same fashion as Equation j3 except with j3fs replaced by 

PcD. 

The pressure in regions 2 and 3 are identical, and are found 
using the strong shock jump conditions. 



P2=P3 = l^PlSM^FS^ 



(14) 



One important aspect to note is that the reverse shock is not 
necessarily strong or relativistic in the frame of the fluid in 
region 4: if the forward shock is relativistic, the reverse shock 
is not (see PW for details). However, if the blast wave en- 
counters a strong density jump causing the forward shock to 
become non-relativistic, the reverse shock can form with rel- 
ativistic speeds. Since we cannot assume a relativistic reverse 
shock, we use the full shock jump conditions yielding the fol- 
lowing equation for the density in region 3: 



P3 = 



4;^373P4 
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(15) 



where 73 is the Lorentz factor of the fluid in region 3 from the 
frame of the unshocked fluid in region 4: 
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Lastly, we need to calculate the velocity and radius of the re- 
verse shock. We use conservation of mass to obtain the veloc- 
ity of the reverse shock: 



PRS = 



P3l3p3-P4l4f^4 



(17) 



P373-P474 

To calculate the radius of the reverse shock, we use conser- 
vation of mass as well. The total mass in regions 3 and 4 is 
constant after the encounter, meaning the integral over these 
regions equals the swept up mass prior to the encounter. 

PRrS PRCD pRenc 

r^p\dr+ r^p',dr= r^p^tUMr (18) 

f Rback ^ RrS ^^ 



Equation [18] is important because we know the densities in 
regions 3 and 4, we have calculated Rback from Equation [8j 



we know Rqd from Equation[T3]but with Pfs replaced by Pcd, 
and we know the total mass swept up prior to the encounter 
(the right hand side of Equation [18]). This leaves us with an 
equation for the radius of the reverse shock: 



Rls(p4-P3)-- 
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^encPwind(Renc) -^ RbackP4~ ^CdP3 ' (1^) 



Equations [Ts] and 19 contain p\ which is the density of the 
fluid in the lab frame, p^ = 7/?, whereas p, without a prime, is 
the co-moving density. 



Solving for Rrs in Equation [19] is only applicable when the 
reverse shock has not yet passe3^he back of the shock. Once 
the reverse shock passes the back of the shock. Equation [18 
no longer applies, and the new integral to solve is: 
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where we now solve for Rback instead of Rrs, and Rback is no 
longer calculated from Equation [8] because region 4 has been 
shocked by the reverse shock ana has knowledge of the en- 
counter. The blast wave continues on with these two regions 
of Equation [2Q|-one containing the newly swept up mass in 
region 2, anoone containing the mass swept up prior to the 
encounter-resulting in the blast wave never completely re- 
turning to the BM solution of a bla st wa ve in an ISM. 

In the following section. Section 2.3 we compare our ana- 
lytical solution with the numerical results for various types of 
encounters. 



2.3. Numerical Results 

Our numerical simulations demonstrate the resulting dy- 
namics of an afterglow blast wave in a stellar wind environ- 
ment encountering a sudden change in density and an ISM 
environment. We use the relativistic adaptive mesh (RAM) 
parallel RHD code ( .Zhang & MacFadyen. 2006, 2009) with 
a second-order weighted scheme. We simulate two differ- 
ent size jumps ("walls"), of factors 4 and 100, and one drop, 
of factor Y^ . We use these specific initial conditions to nu- 
merically simulate a wide range of encounters thought to po- 
tentially cause light curve flares. However, as we explain in 
Section [3] none of these scenarios result in flares on the time 
scale of 0.1 to 10 days. Before delving into the resulting light 
curves, we first explain the dynamics and numerical results. 

We set up six simulations all of which have a starting time 
and radius equivalent to the fluid Lorentz factor behind the 
front of the shock, 7 = 15. This starting Lorentz factor was 
chosen to ensure 7 > ^, where Oo = 0.l rad (typical for after- 
glow jets, ^axT^^aLpOOT]) for our simulations and is the half 
opening angle of the jet. 
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We first ran all six simulations in one dimension, and then 
ran three simulations in two dimensions, choosing to run those 
which would be most informative to answering the questions 
outlined in the Introduction. 

Figure [2] shows the analytical solution of the density of sim- 
ulation 75a4, in one dimension plotted against the numerical 
result of simulation 75a4 in one dimension, showing the cor- 
respondence of our analytical model with our numerical re- 
sults of the shock formation at the encounter. Our analyti- 
cal model assumes homogeneous shells with constant density 
between the shocks which is why the analytical solution in 
Figure [2] has a clearly marked "back" of the shock, whereas 
the numerical simulation uses the full BM profile and has no 
exact "back" of the shock. The small discrepancy between 
the analytical and numerical Rfs shown in Figure |2]is roughly 
ARfs = 9.6 X 10^^ cm. Our analytical solution is accurate 
to the order of 1/7^ = 0.0855 at the time shown in Figure [2] 
and ^^ = 3.893 X 10"4 < 1/7^ meaning that the analytical 
model conforms to the numerical results to the required ac- 
curacy. Also, it is apparent from Figure [2] that the densities 
within each region are not constant, but actually have a posi- 
tive slope. This is a known discrepancy between the assumed 
homogeneous shells and the full BM solution, and leads to a 
slight over-predictio n of the flux during th e encounter by the 
analytical solution ( Van Eerten et al.|2009| ). 

Figures [3] and [4] show our analytical solutions plotted with 
the simulation results of the Lorentz factors behind the for- 
ward shocks as they evolve. The small deviations towards the 
beginning of the simulations (at smaller radii) in these figures 
are from lack of resolution. The resolution required to resolve 
the blast wave at these early times is very high. However, for 
the simulations with encounters at 72 = 5, the discrepancy dis- 
appears well before the encounter. For the simulations with 
the encounter at 72 = 10, the lack of resolution does affect 
the convergence of the numerical and analytical results, but 
qualitatively, the results are the same. The resolutions of our 

The devia- 
are 



simulations are discussed further in Appendix \K 

tions at large radii (at very late time) in Figures 

from the Lorentz factor of the forward shock becoming"non 

relativistic, and our analytical model having accuracy of order 
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Fig. 2. — Figure of the shocks formed during the encounter for simulation 
75fl4 in ID (wind termination shock at jenc = 5, and a density increase of 
factor 4). The simulation is plotted in blue, the dashed black line shows the 
analytical solution for this simulation, and the green dashed line represents 
simulation no encounter in ID (no circumburst environment change) at the 
same time as simulation 75^4 is depicted. 
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Fig. 3. — Lorentz factor of the fluid behind the forward shock plotted over 
radius for the analytical solution as well as the numerical simulation for sim- 
ulation 75fl4 (wind termination shock at ^enc = 5, and a density increase of 
factor 4) and 710^4 (wind termination shock at jenc = 10, and a density in- 
crease of factor 4) in ID. 



From the examples shown in Figures [2] through |4j it is ap- 
parent that our analytical solution accurately predicts the dy- 
namics of a one dimensional blast wave in a stellar wind en- 
countering a change in density, be it a jump or drop, followed 
by an ISM environment. Also, it is important to note that 
when the post encounter forward shock remains relativistic, 
the fluid Lorentz factor at which the encounter occurs does 
not qualitatively change the post-encounter dynamics in one 
dimension, as seen in Figure |3] The post-encounter dynamics 
are simply scaled to a reflect the larger encounter Lorentz fac- 
tor. This is also clear from the analytical model, since neither 
the change in Lorentz factor at the encounter nor the pre- and 
post-encounter slopes depend on Lorentz factor. 

The one dimensional simulations show that the dynamics 
follow the analytical predictions that the Lorentz factor of the 
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Fig. 4. — Lorentz factor of the fluid behind the forward shock plotted over 
radius for the analytical solution and the numerical simulation for simulation 
75flT^ in ID (density decrease of factor 100 at jenc = 5). 



forward shock drops immediately as it encounters a higher 
density region, and the amount by which it drops is propor- 
tional to the change in density at the encounter. This is shown 
in Figure0for simulations j5a4, j5a-^, and 7lOalOO. Sim- 
ulations 75^100 and jl0a4 show similar behavior. The Left 
Column in Figure [5] shows the Lorentz factor over time and 
radius for simulations j5a4 and 7 10a 100 in ID. If the blast 
wave encounters a drop in density, the Lorentz factor of the 
blast wave will increase (i.e., the blast wave will speed up), 
again, in proportion to the size of the drop, which is seen in 
the Top Right Figure in Figure [5] 

The effects of the magnitude of density change on the dy- 
namics of the reverse shock are shown in Figure [5] The Top 
Left Figure in FigurelS] depicts the Lorentz factor for the sim- 
ulation with a small density increase of factor 4. Shortly after 
the encounter, the reverse shock slowly travels away from the 
forward shock towards the origin, but this is not an immediate 
process like that of a larger density jump shown in the Bot- 
tom Left Figure in Figure [5] For a large density jump, almost 
immediately after the encounter, the reverse shock has a large 
enough Lorentz factor, compared to the Lorentz factor of the 
forward shock, to break away from the forward moving fluid 
and travel backwards towards the origin. This difference in 
reverse shock Lorentz factors of simulation 'y5a4 and simu- 
lation 710^100 is from the forward shock becoming immedi- 
ately non-relativistic after the encounter with a large density 
jump (of factor 100) versus the much smaller decrease in for- 
ward shock Lorentz factor with a smaller density jump (of 
factor 4). As explained in Section 2.2 the amount by which 
the forward shock Lorentz factor drops is proportional to the 
Lorentz factor with which the reverse shock forms. For sim- 
ulation 75a Y^, which has a density decrease of factor 100, 



the reverse shock takes much longer to overcome the forward 
moving fluid and start traveling backwards towards to origin. 
This is seen in the Top Right Figure in Figure [5] The Lorentz 
factor of the forward shock increases, causing trie Lorentz fac- 
tor of the reverse shock shortly after the encounter to be very 
small in comparison with the forward shock Lorentz factor. 

A closer look at the Bottom Left Figure in Figure [5] reveals 
the reflecting-like behavior of the inner boundary as the re- 
verse shock travels towards the origin and is seemingly re- 
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fleeted baek. This behavior is independent of boundary condi- 
tions or inner boundary position and results from the increase 
in pressure caused by the reverse shock traveling towards the 
origin. This increase in pressure creates conditions similar 
to a fireball, causing a new blast wave to be formed. This is 
not seen with simulations 75^4, 710^4, or J^a^ which have 
smaller density jump factors/a density drop, but is seen with 
simulation 75^100 (wind termination shock when 72 = 5 with 
a density jump of 100) yielding the conclusion that even at 
slower speeds, the strong, large density jump still slows the 
forward shock enough to cause a strong reverse shock to form 
and in turn, cause the conditions similar to a fireball. From 
here forward in this paper, we will refer to this newly formed 
blast wave caused by the strong reverse shock as the "sec- 
ondary blast wave" for simplicity and clarity. 

The formation of the secondary blast wave is also seen in 
our two dimensional model of simulation 7 10a 100. How- 
ever, the fluid instabilities that follow the encounter with the 
wind termination shock in two dimensions slow the secondary 
blast wave causing it to never encounter the wind termination 
shock, and in turn, the continuous reflecting back and forth 



is never realized. This is shown in the Bottom Right Figure 
in Figure |5] which displays the fluid Lorentz factor over time 
and radius of a one dimensional slice of the two dimensional 
simulation 7 10a 100 along the axis of the jet. 

In order to study the effect of the density change on the post 
encounter blast wave collimation, we ran simulations 75a4, 
7 10a 100, and 75a j^ in two dimensions. Figure pi shows the 
results of our two dimensional simulations. We fina that when 
the blast wave encounters a sudden density increase that is 
not very large (of factor 4, like that of simulations 75a4 and 
7l0a4), the fastest, most energetic fluid will puncture the wall 
and continue through to the new medium (Top Left Figure in 
Figure [6]). This behavior results in the blast wave recollimat- 
ing after the encounter as only the fluid with the highest pres- 
sure can pass through to the new region and open a pathway 
for the rest of the blast wave to continue through to the new 
environment. This is most clearly depicted in the Top Right 
Figure in Figure |6] which shows the percentage of energy en- 
closed at various angles over time for simulation 75a4 in two 
dimensions. Before the encounter, the blast wave spreads out, 
which is seen as the angles encompassing certain percentages 



of energy increase over time at times less then Tenc in the 
Top Right Figure in Figure [6] At the time of the encounter, 
the most energetic fluid passes through to the new environ- 
ment and recolHmates, which compresses a large percentage 
of the blast wave energy into a smaller angle. This accounts 
for the decrease in angle encompassing roughly half of the 
total energy after the encounter in the Top Right Figure in 
Figure |6] Also, if there were no encounter, the angle encom- 
passing roughly 80% of the total energy would continue to 
increase as it did prior to the encounter in the Top Right Fig- 
ure in Figure[6j but instead, the angle stays relatively constant. 
This is because the fluid of higher energy (depicted with col- 
ors of yellow through red in the Top Left Figure in Figure |6]) 
no longer spreads sideways but tries to travel through the wall 
to the new medium. The only part of the blast wave that con- 
tinues to spread sideways is the lowest energy fluid (shown 
in green and light blue in the Top Left Figure in Figure [6]), 
which results in the continual increase of the angle encom- 
passing 100% of the total energy of the blast wave for simula- 
tion 75a4 in 2D (shown in the Top Right Figure in Figure |6]) 
after the encounter. 

The dynamics of a blast wave encountering a large density 
increase (on the order of 100, like that of simulations 75^100 
and 7 10(2 100) are different than that of a smaller density in- 
crease. The most energetic fluid punctures the wall and con- 
tinues on to the next medium, similar to the most energetic 
fluid encountering a smaller density increase, but instead of 
opening up a pathway for the rest of the blast wave to move 
through, a strong reverse shock is formed. This strong reverse 
shock shocks and spreads the fluid that has not yet passed 
through to the new medium, forming vortices that are trapped 
within the stellar wind environment and never pass through 
to the ISM. The Middle Left Figure in Figure [6] and Fig- 
ure |7] illustrate these dynamics showing the specific energy 
of simulation 710^100 in two dimensions at two times after 
the encounter. The Middle Left Figure in Figure [^depicts the 
blast wave a short time after the encounter showing that the 
most energetic fluid has passed through to the ISM, and the 
reverse shock is traveling back through the fluid in the stel- 
lar wind, shocking and spreading the fluid. At the much later 
time of Figure 12 the highest energy fluid that passed through 
to the ISM hasbegan to expand as it assimilates to the new 
medium. Also, the reverse shock has now passed through the 
fluid behind the density jump, caused the formation of vor- 
tices, and created the condition for the secondary blast wave 
to form. The remains of secondary blast wave can be seen 
in Figure [7] by the vorticies near the horizontal axis at radii 
around 2 x 10^^ cm. 

The Middle Right Figure in Figure [6] shows the percentages 
of energy encompassed in various angles, depicting that the 
highest concentration of energy is contained within the angle 
holding the fluid that passes through the density jump, and a 
very minimal amount of energy is contained within the vor- 
tices that have been formed in the stellar wind region. The 
chaotic behavior of the region in purple and pink in the Mid- 
dle Right Figure in Figure|6]is from those vortices as well. 

Lastly, the dynamics of tne simulation with the sudden den- 
sity drop do not show a recollimation with the increase in 
Lorentz factor. The fluid speeds up and as it does, it spreads 
at a rate slightly faster than in the stellar wind environment. 
This is shown in the Bottom Row of Figure [6] 

With the understanding of the dynamics of the blast wave 
encountering a sudden change in density, the next section dis- 



cusses the resulting light curves. 

3. LIGHT CURVES 

We calculate the GRB afterglow ligh t curves using the 
radiation calculation methodology from [Van Eerten et al.| 
( 2010a|b ) at an X-ray frequency of 5 x 10^ '' Hz (except where 
otherwise specified), similar to the frequency detected by 
Swift. The time at which the flux from the local emission of 
the afterglow is observed is denoted as the "observer time," or 
tobs^ and is calculated by: 



''obs ■ 



rcosO 



(21) 



where te is the emission time (i.e., the lab frame time of the 
emitting fluid element), r is the radius of the local fluid ele- 
ment at te, and 6 is the angle between the direction to the fluid 
element and the direction to the observer. For an emission 
time, te, the fluid at the front of the shock will be observed 
first, and can be represented as: 



-^obs — te 



R 

c 



(22) 



where Tobs denotes the earliest time any flux from the blast 
wave at the time, te, is observed. On our figures, we label the 
first time the encounter can be observed as Tobs,enc- 

For our light curves, we use values of eg = 0.01, e^ = 0.1, 
and p = 2.5, where eg and e^ are the fractions of internal en- 
ergy that contribute to the magnetic field at the shock front 
and to accelerating electrons respectively, and p expresses the 
slope of the energy distribution of the shock- accelerated par- 
ticles (and is not to be confused with pressure). We also do 
not consider electron cooling in our light curves shown here, 
but have found that this does not qualitatively change our re- 
sults. We expect the pre-jet break light curves to follow the 
temporal behavior of ^^^"^^^/"^ in the stellar wind environment, 
and r^(i-/^)/4 in the ISM (e.g., Granot & Sari 2002 ). 

It has been claimed that blast waves encountering wind ter- 
mination shocks generate flares or rebrightenings which can 
be seen at early times and are caused from either the reverse 
shock forming at the encounte r (PW) or the b last wav e transi- 
tionin g to the new medium (E ldridge et al.l20 06i Mesler et al. 
2012| ). As discussed below, from our analytical and numeri- 
cal analysis,we find that a GRB jet encountering a change in 
circumburst medium does not cause a flare of the kind seen 
by Swift. Neither the reverse shock at the encounter, nor the 
blast wave's transition to the new medium at the encounter 
will cause an observable flare. 

3.1. No Flares During the Observed Time Scale 

Figure [8] depicts that there are no observable flares in the 
flux emitted from a blast wave encountering a small jump, 
large jump, or large drop. This figure shows the light curves 
calculated from the one dimensional simulation and analytical 
model as well as the light curves calculated from the two di- 
mensional numerical simulations for on and off-axis observer 
angles. For all cases studied in this paper, there were no ob- 
servable flares. 

The simulations demonstrate that a blast wave encountering 
a wind termination shock of a small increase in density (of 
factor 4) does not result in any disruption in the flux observed 
at early times. In fact, for a jump of factor 4 (simulations 
j5a4 and 7 10^4), the light curve gradually transitions to the 
slope expected in the new environment, confirming the results 
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from the spherical cases of NG and | Van Eerten et aL]p009| ). 
This is seen with the on-axis light curves of the one and two 
dimensional simulations (Top Left in Figure [8]) and with the 
light curves for a wide range of observer angles of the two 
dimensional simulation 'y5a4 (Top Right in Figure [8]). The 
encounter is first observable at roughly 46 days for the simu- 
lations with jenc = 5, and for our simulations of jenc = 10, the 
first time the encounter is observed, Tobs.enc ^ 3 days. These 
times are much later than any flares observed from Swift and 
follow from our choice for 7^^^, which was selected for nu- 
merical reasons. However, the time at which the encounter is 
observed is inversely proportional to V\ 



" enc 



^obs.e. 



(23) 



This scaling is calculated using Equations |6] and |7j and noting 
that at (and before) the encounter k = 2. It follows that mod- 
eling encounters at larger Lorentz factors, of say '^enc = 25, 
which corresponds to Tobs,enc ^ 0-07 days, would result in 
the observed encou nter being in the ran ge of times flares are 
observed (see e.g., jZhang et al.||2006| ). Nevertheless, from 
the Top Right Figure in Figure M it is apparent that the light 
curves for blast waves encountering changes in external den- 
sity profiles at varying Lorentz factors scale for one dimen- 
sional dynamics. Thus the light curve behavior at the time of 
the encounter for our numerical simulations of jenc = 5 and 10 
are indicative of the behavior of encounters at larger Lorentz 
factors. To ensure that this is accurate, we analytically model 
encounters at jenc = 25, and discuss the light curves from those 
dynamics throughout this section. The light curves calculated 
from our analytical model use the same linear radiative trans- 
fer algorithm as the light curves calculated from the numerical 
simulations. 

For simulation 75^2^^, the flux observed follows the same 
qualitative results as simulation 75a4 at early times-the ob- 
served flux transitions smoothly to the slope of the flux ob- 
served solely in the ISM environment after the encounter. 
Even though the blast wave's velocity increases after the en- 
counter (Figures]?] and Top Right in|5]), a flare during the en- 
counter is not observed. Instead, the flux drops more steeply 



as it evolves towards its lower post-encounter base level. This 
is shown in the two figures in the Bottom Row in Figure [8] 
where the figure on the Bottom Left depicts the flux observea 
on axis with the blast wave using the one and two dimensional 
simulations along with the analytical solution, and the latter 
displays the flux observed at various observer angles for the 
two dimensional simulation. 

It is apparent from the light curves of simulations 75a4 and 
75a Y^ , that as the magnitude of the flux of a blast wave in the 
ISM environment lowers, there is no rebrightening. Next, we 
discuss the light curves of simulation 7 10a 100, which are of a 
blast wave encountering an ISM environment of much higher 
density. 

With the magnitude of the flux of a blast wave traveling 
solely in the ISM environment being much higher than the 
magnitude of the flux from a blast wave in the wind (because 
the magnitude is a function of the density in that region), one 
might first assume a flare must be observed as the flux from 
the blast wave at the encounter evolves to the new slope and 
magnitude. However, this is not correct-there is still no re- 
brightening at the time of the encounter. This is shown in 
the Middle Row of Figure [8] The off axis light curve calcu- 
lated at Oobs = 0.2 radians m the Right Figure in the Middle 
Row of Figure [8] shows a bump at around 100 days. Although 
the flux does not increase at this point, it does deviate from 
the expected slope. This could be considered a rebrightening, 
however, it is only seen at observer angles larger than the jet 
opening angle, and it is unlikely that the prompt emission will 
be observed at this angle due to the extreme beaming needed 
of the gamma radiation. 

Figure [9] depicts the transition of the light curve as the blast 
wave changes environments as well as the overall scaling of 
the light curves. The Top Left Figure in Figure |9] shows ana- 
lytical light curves for a blast wave traveling in a stellar wind 
environment with '^enc = 25 for two scenarios: a density jump 
of factor 100 (the solid blue line) and a density drop of factor 
10000 (the dashed blue line). These light curves are plotted 
against the light curve for a blast wave traveling solely in the 
stellar wind environment (green line), a blast wave traveling 
solely in the homogeneous medium corresponding to the den- 
sity jump of factor 100 (red solid line), and a blast wave trav- 
eling solely in the homogeneous medium corresponding to the 
density drop of factor 10000 (dashed red line). As seen in this 
figure, the observed flux after the encounter slowly transitions 
to the slope of the flux observed from the ISM environment, 
resulting in no sudden increase in flux. The flux does not im- 
mediately jump up to the new magnitude because the blast 
wave has been slowed by the encounter and is not energetic 
enough to cause a rebrightening. For the case of a drop-the 
flux observed from a blast wave in the new ISM environment 
is at a lower magnitude than the flux observed from a blast 
wave in the wind environment. This results in the flux ob- 
served from the blast wave encountering the new environment 
to decrease, instead of increase, after the encounter to evolve 
in the new environment. 

In addition, this figure shows an interesting feature, that can 
also be seen in Figure [8j of a shallowing of the light curve. 
The light curve transitions from the slope of a blast wave in 
a stellar wind to the less steep slope of a blast wave in an 
ISM. This shallowing of the light curve is an observed feature 
of GRB afterglows, and a change in circumburst medium is 
one of the few ways for the light curve of the blast wave to 
transition to a new slope. 
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From the light curves of simulations j5a4, 7 10a 100, and 
75a j^, as well as the light curves from the simulations not 
shown here (75^100 and 7 10^4), we conclude that there is 
no rebrightening at times on the order of minutes to days for 
a blast wave encountering a change in circumburst medium. 
The simulations show that when the magnitude of the flux ob- 
served in the ISM is lower than the flux observed in the stellar 
wind environment at the time of the encounter (Top Left Fig- 
ure in Figure |9]) the blast wave conforms to the new slope and 
there is no rebrightening. When the magnitude of the flux 
observed from a blast wave in an ISM environment is higher 
than the magnitude of the flux observed from a blast wave 
in the stellar wind environment, the blast wave still smoothly 
transitions to the new environment. In all cases, though, the 
flux observed from the blast wave after the encounter is lower 
than the flux observed from a blast wave that has only traveled 
in the ISM. 

To test whether the light curves for our simulations are an 
accurate depiction of a blast wave encountering a circumburst 
density change, or if they only represent the light curves for 
encounters at fluid Lorentz factors of 7 < 10, with changes 
in density of factors < 100 and with 6>o = 0.1, we analytically 
solve for the light curves of various blast wave encounters for 
various encounter Lorentz factors, jenc^ density change fac- 
tors, and jet half opening angles, ^o- The Top Right Figure 
in Figure |9] shows light curves for one dimensional analytical 
models ofblast waves traveling in a stellar wind environment 
that encounter a density change of factor 10 followed by a ho- 
mogeneous medium with the encounters happening at vary- 
ing Lorentz factors. When the factor by which the density 
changes is kept fixed, the shape of the light curve remains in- 
dependent of the actual time of the encounter and light curves 
for different encounter times can simply be shifted to match, 
as can be seen from the figure. 

The Bottom Left Figure in Figure |9] shows the light curves 
for blast waves traveling in the stellar wind environment that 
encounter a density change of varying magnitudes at jenc = 25. 
We show here density changes ranging fr om 10~^ t o roughly 
correspond with the density drops used in Mesler et al. (2012), 
to 10^. For larger density jumps, the blast wave quickly be- 
comes non-relativistic, causing our analytical model to be- 
come less accurate. Lastly, the Bottom Right Figure in Fig- 
ure [9] shows light curves for blast waves traveling in a stel- 
lar wind environment that encounter a density jump of factor 
10 at jenc = 25 for various jet half opening angles, ^o- The 
only light curve in this figure that has some resemblance of a 
flare is that of the blue line, which has an opening angle very 
close to jenc ^ V^o, meaning that for any flare to possibly oc- 
cur around the time of the encounter, jenc = l/^o- This is the 
setup of our simulations 7 10^4 and 7 10a 100: jenc = 10, and 
^0 = 0.1, and no flare is observed in one or two dimensions. 

From this analysis, we conclude that for a blast wave trav- 
eling in a stellar wind environment that encounters a sudden 
change in circumburst environment, there is no observable re- 
brightening regardless of the size of the jump, drop, or fluid 
Lorentz factor at the time of the encounter that we have stud- 
ied. There is a bump in the light curve at the observation angle 
of Oobs = 0.2 radians for the light curve of Simulation 7 10a 100 
in 2D, but Oobs > ^0 and thus this is not a likely indication of 
an observable flare. We do not see a flare from the forward 
shock, nor do we see a flare caused by the reverse shock at 
times corresponding to the afterglow flares seen by Swift. 

We study a single spherically symmetric change in the cir- 



cumburst environment. However, the actual environment may 
have a more complicated structure with instabilities (e.g.. 



van 



Marie & Keppens]|2012) or thin circumburst shells ( Mesler 



et al.||2012| ). A study of asymmetric structures (i.e., density 
clumps) around the progenitor would require three dimen- 
sional simulations and is beyond the scope of this paper. We 
do look at the interaction with circumburst shells, and our re- 
sults are discussed further in Appendix [B] In this case too, we 
have been u nable to reproduce the flares reported by [Mesler 
[eFaLldMni ). 

3.2. Late Time Rebrightening 

The simulations of a density jump did show flares at times 
long after the encounter occurred (at times much later than 
plotted in Figure [s]). At such late times (on the order of 10^ 
days), these flares are beyond X-ray observations, but may be 
of interest at radio frequencies. There are two contributors to 
these flares-the fireball conditions created from the strong re- 
verse shock that causes a secondary blast wave to form, shown 
in the Bottom Row of Figure [5] and the encounter of the re- 
ceding jet with the wind termination shock. We discuss the 
contribution of the secondary blast wave first, followed by a 
discussion of the receding jet. 

Figure [To] shows the very late time flare observed from the 
one and two dimensional simulations of simulation 7 10a 100. 
Comparing the time at which the flare begins for simulation 
7 10a 100 in ID and 2D, shown in Figure 10 with the time 
and radius of the secondary blast wave shown in the Bottom 
Left Figure in Figure |5] it is apparent that the secondary blast 
wave is a cause of the rebrightening. This is confirmed by 
the Left Figure in Figure [T0| which shows the light curve cal- 
culated for simulation 7 10a 100 in ID using the entire blast 
wave, against the light curve observed from the part of the 
blast wave that is at a radius larger than the encounter radius. 
The Left Figure in Figure [TO] clearly illustrates that the flare is 
caused by the fluid insideme encounter radius at late times, 
and the most energetic fluid in the ID simulation that is at 
radius smaller than Rene at late times is the secondary blast 
wave. 

The Left Figure in Figure [T0| although explaining the flare 
for the one dimensional case of simulation 7 10a 100, does not 
explain the small flares seen in the one dimensional conical 
outflow for simulation 75a4 nor the flares for the two dimen- 
sional simulations 75a4 and 7 10a 100. 

To fully understand the cause of this late time flare, we must 
emphasize that it is seen at a wide range of observer angles 
and frequencies (Right Figure in Figure [TO]) at the same ob- 
server time, and is observable at the same time as the flare 
from the secondary blast wave. The secondary blast wave 
cannot be the cause of the flare for the 2D simulation because 
it is quickly slowed by the surrounding material (seen in the 
Bottom Right Figure in Figure [5]), and is also not a character- 
istic of simulation 75a4 in one or two dimensions. Moreover, 
the flare observed in simulation 75a4 is observable in the 2D 
light curves and the ID light curve that assumed conical out- 
flow, but not in the ID light curve assuming spherical outflow. 

The cause of the late time flare is the receding jet that en- 
counters the sudden density jump (but not from a density 
drop). The reason this event seemingly occur at the same 
time as the secondary blast wave is because the sum of the 
distance traveled by the reverse shock and the flux from the 
secondary blast wave is the same as the distance traveled by 
the flux observed from the receding jet encountering the den- 
sity jump. In short, the density jump results in a sudden in- 
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X-Ray Light Curves, Simulation 710a 100 in ID 



Light Curves, Simulation 7lOalOO in 2D 



— All Jet included 

— Including Jet>Rc] 




10^ 10^ 10^ 

Observer Time (days) 



10^ 




10^ 10^ 10^ 

Observer Time (days) 



Fig. 

of6>n 



10. — Figures showing the late time flare of simulation 710^100 (^yenc = 10, observed at roughly 3 days, with a jump of 100 and the jet half-opening angle 
0.1). Left: light curve of the one dimensional simulation 7 10^100 including the entire evolution (blue line) and including only the fluid at radii larger 



than Rene (green line). Right: light curves for the 2D simulation 710^100 at three frequencies: radio, u = 4.86 x 10^ Hz, visible, u = 4.56 x 10^^ Hz, and X-ray, 



: 5 X 10^^ Hz. These light curves were calculated on axis with the blast wave, 
radio frequency using only the forward jet (and not including the receding jet). 

crease in flux traveling in the opposite direction of the blast 
wave itself, meaning a spike in the light curve will not be ob- 
served at the time of the encounter because the light from the 
encounter is traveling away from the observer. This is con- 
firmed by the green lines in the Right Figure of Figure [T0| 
which shows that the light curve observed only from the for- 
ward jet does not result in a very late time flare whereas the 
light curve observed from both the forward and the receding 
jet does have a flare. We obtain the same results when an- 
alyzing simulation 7 10a 100 at different frequencies, and by 
analyzing the very late time flare of simulation 'y5a4 in 2D 
(and ID conical outflow). 

From our simulations, the only flare calculated is that at 
times much later than flares observed by Swift and is at a 
much lower magnitude in flux than Swift observes. The only 
possible flare on the time scale of observed flares we are able 
to observe with a sudden change in circumburst medium is 
the small bump in the light curve seen at around 30 days at 
Oobs = 0-2 radians in the Middle Right Figure in Figure [8| but 
this occurs only for Oots > ^0 and thus is unlikely to be ob- 
served. 

4. SUMMARY AND CONCLUSIONS 

We have shown numerically and analytically that a blast 
wave evolving partially in a stellar wind environment that en- 
counters a sudden change in density, either an increase or a 
decrease, followed by a constant density environment for a 
wide range of initial conditions does not cause an observable 
rebrightening. Flares at very late times are a function of the 
size of the density jump-the larger the density jump at the en- 
counter, the stronger the reverse shock, resulting in a brighter 
late time flare. We studied light curves for a wide range of 
frequencies, and although this paper focuses on X-ray light 
curves, we saw similar behaviors for a wide range of frequen- 
cies. These flares potentially might be observable at radio 
frequencies under favorable conditions (Figure 10). This an- 
swers the first question listed in Section fl} the size of the den- 
sity jump does affect the dynamics andresulting light curve 
and may cause a faint very late time flare at radio frequencies 
for large jumps, but there are no observable flares at X-ray 
frequencies. 



^obs = 0.0 radians. The blue dashed line shows the light curve calculated at a 

We found that for a blast wave traveling in a stellar wind 
environment encountering an ISM environment, the resulting 
flux observed will gradually transition from one environment 
to the next. As seen in Figure |9] if the flux observed from a 
blast wave traveling solely in tne ISM is lower than the flux 
observed with a blast wave traveling solely in the wind envi- 
ronment at the time of the encounter, the flux observed from 
the blast wave will simply dim and follow the same light curve 
slope of the ISM. If the flux observed from a blast wave trav- 
eling solely in the ISM is higher than the flux observed from a 
blast wave traveling solely in a wind environment at the time 
of the encounter, the observed flux from the blast wave will 
not suddenly increase, but stay relatively steady as it transi- 
tions to the new slope of the ISM. This analysis was mostly 
done for the one dimensional simulations and showed that the 
light curve is affected by the size of the density change, but 
not the Lorentz factor at the encounter. We have explored var- 
ious ratios of'^enc/ljet with two dimensional simulations and 
found that this conclusion holds in two dimensions. How- 
ever, strictly speaking the recollimation seen in the two di- 
mensional simulations of density jumps is dependent on the 
Lorentz factor of the fluid at the encounter, and the amount 
by which the blast wave recollimates should be inversely pro- 
portional to the fluid encounter Lorentz factor (i.e., blast wave 
recollimates more for lower fluid encounter Lorentz factors), 
although this is not expected to have an effect on the occur- 
rence of a flare. 

We have studied the two dimensional and one dimensional 
effects of a density drop and have shown that the blast wave 
does increase in speed, but does not recollimate. There is also 
no rebrightening caused by this sudden increase in blast wave 
speed (Figure!^. Our two dimensional studies of a density 
jump have yielaed the conclusion that there is some sideways 
spreading from a high energy collision of the blast wave with 
a large jump in circumburst density, and the amount by which 
the blast wave spreads is highly dependent on the size of the 
density jump. There is also a flare at very late times, but this 
is not caused by the sideways spreading. The late time flare is 
caused from the reverse shock of the receding jet (Figure 10). 
This very late time flare only occurs years after the initial flux 
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observed from the blast wave, and data from Swift is only 
gathered for days after the first flux is observed resulting in 
the conclusion that this is not the flare observed by Swift. In 
addition, these flares range from being 10~^ to 10~^^ orders of 
magnitude smaller than the initial flux of the blast wave which 
is below the threshold of Swift. 

We have analytically modeled the light curves of one di- 
mensional simulations, and are able to use the model for on 
axis two dimensional dynamics that do not cause strong re- 
verse shocks as well. From our analytical solutions and our 
numerical models, we have answered the questions listed in 
the Introduction and have concluded that a blast wave travel- 
ing in a stellar wind environment that encounters a change in 
density followed by an ISM environment will not cause ob- 
servable flares. 

Our work confirms earlier studies (NG, |Van Eerten et al.| 
[2009 ) that a wind termination shock will not cause flares in 
light curves and extends this analysis to density drops and 
to two dimensional simulations. We were unable to repro- 
duce the flares seen in PW due to the discrepancy of the fluid 
Lorentz factor during the encounter-the analysis done in PW 
assumed a constant fluid Lorentz factor during the encounter, 
and we have extended this analysis to more accurately model 
the decrease in fluid Lorentz factor at the encounter. We were 



not able to reproduce the flares seen in [Mesler et al.| ( |2Q12| ) 
and |Eldridge et al.| ( |20Q6| ). Both of these studies use an ana- 
lytical model for the blast wave based on total swept up mass. 
They do not account for the radial or angular structures of 
the blast wave, such as the differing shock regions shown in 
Figure [T] and the spreading of the shock as it becomes non- 
relativistic. Without accounting for these radial and angular 
structures, their models of the change in Lorentz factor with 
change in mass result in a discrepancy in time scale and flux 
observed between ours and their results. In addition, as the 
blast wave becomes non relativistic, the blast wave spreads 
outwards and its width becomes comparable to the radius of 
the blast wave. Thus, for the the interaction of a blast wave 



with a spike, like that in Mesler et al. (2012), the swept up 
mass model does not allow for capturing features on the same 
or shorter time scale of the blast wave crossing the shell be- 
cause the width of the shell is smaller than the width of the 
blast wave. This is discussed more thoroughly in Appendix[B] 
The very late time radio flares from our simulations are sim- 



ilar to those seen in s imulations of slower ejecta catching up 
with the blast wave ( [Maxham & Zhang||2Q09[ |Vlasis et"aL 
|2011| ). However, the simulations of the slower ejecta catch- 
mg up with the blast wave are simulated by adding energy to 
the system at late times leading to a different light curve slope 
after the encounter. The very late time radio flares seen in 
this paper were calculated from a simulation in which no en- 
ergy was added, resulting in the flares seen in this paper being 
of different origin than from slower ejecta reaching the blast 
wave. 

We conclude that a wind termination shock, or more gen- 
erally, any sudden transition in circumburst density (even ex- 
treme changes), is very unlikely to be the cause of the flares 
observed by Swift because those flares r eturn to the same base - 
line as the light curve prior to the flare ( [Burrows et al.|2005] ). 
With the encounter of a change in circumburst environment, 
the only observable flares do not return to the same base- 
line because flux from the receding jet is also being observed. 
Moreover, if the flare occurred at the time of the encounter, the 
light curve would transition to the new slope of the ISM mate- 
rial, not the slope of the stellar wind. This makes it highly un- 
likely that a flare seen at the encounter with the change in cir- 
cumburst environment could be the explanation for the flares 
seen by Swift. 

Through our numerical and analytical analysis of a blast 
wave encountering a circumburst density change that origi- 
nated from the same point, we have concluded that a flare will 
not be observed by Swift, but a change in environment could 
be the explanation for the shallowing observed for some GRB 
afterglows ( Evans et al.,20Q9) . It would be of interest to com- 
pare our results with observational data in future research. 
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APPENDIX 
GRID SETUP AND RESOLUTION 

To run our simulations, we first found the radius and time corresponding to the Lorentz factor at which we wanted the external 
density change to occur, the size of the grid, as well as the needed resolution of the simulation. We use Equations pi and [5| to set 
up the initial shock location and time respectively with the initial fluid Lorentz factor of 15. Next, we use Equation|6]to find the 
location of the Lorentz factor at which the jump occurs, and a form of Equation [7] to calculate the time at which the Lorentz factor 
is equal to that at the jump. 

The grid size needs to be large enough to ensure the simulation would run a satisfactory amount of time after the encounter. 
We choose a maximum observation time, tobs^ for the simulation to run, which gives an equation for the maximum radius and 
simulation time. 

D 

tobs=t-- (Al) 



R in Equation Al is the maximum size of the grid, which we will refer to as Vmax^ and can be represented as the sum of the 
radius of the encounter, Rene, plus some distance covered, Ar. 



R-. 



■ Rene + Ar 



(A2) 



The variable Ar in Equation A2 can be found by integrating the velocity over time. We use a much simplified model for the grid 
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Evolution of 7 for Simulation 75o4 in ID, 
with varying resolutions 




Time (s) 



Fig. 11. — Plot of the Lorentz factor over time for simulation 75^4 (blast wave encountering a density jump of factor 4 at a time and radius corresponding to 
the front of the shock fluid Lorentz factor, ^enc = 5) in ID for various resolutions. 



setup because we need only a rough approximation. Using this simplified model, we assume that immediately after the encounter, 
the blast wave returns to a BM self similar solution of a blast wave in a homogeneous environment. This, again, is not an accurate 
assumption, but accurate enough for this purpose. 



Velocity is represented by v = c/3 and /^ = \/ 1 ~ f- where the subscript of 2 here denotes that this is the Lorentz factor after the 
encounter. T2 if found by a similar equation to Equation [7] except with ^ = 0. 



T2 = Te 



Using the above equations, we obtain the integral of velocity over time. 



/:; 



1 



dt. 



To solve Equation |A4j we can expand the integrand 

3' 



1- 



1 

^ enc 



^enc 



= 1- 



^^ enc 



^enc 



8r4 

^ enc 



(A3) 



(A4) 



(A5) 



Using Equation A5 in Equation A4 to find Ar, we derive Ymax from Equation A2 Next, we calculate the corresponding 



simulation time at the end of the grid using Equation [AT 



We set the minimum value of the grid to a size that ensures the entire initial shock is shown, r^/„ = O.OIT^q. 

Lastly, to find the optimal resolution, we need to be certain that the initial shock is resolved. We use adaptive mesh refinement 
and we specify the num ber o f base blocks, bb^, the number of cells per block, cpb, and the maximum refinement level, ref, for 
the simulation. Equation |A6| describes the relationship these values have to the resolution where A is the amount ISR is refined. 



'max 'm 



AR 



(A6) 



bb;, X cpb X 2ref-i 

AR = Y^, cpb = 8, and the optimal value of A can be tested, as shown in Figureflll Figure fill shows that even with a very small 

amount of resolution in the radial direction, the simulation still converges at the time of the encounter-the part of the simulation 
on which we are focusing. The refinement level in the direction must also be computed. The simulation only computes from 
0<6>< f, which yields: 

^/2 ^ ^0 

bb^ X cpb X 2^ef-i B 

where bb^i is the number of base blocks in the direction, ^0 is the jet half opening angle, and B is the number by which we want 
the width of the jet opening to be refined. 



(A7) 
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Fig. 12. — Figures for the setup with a spike in the circumburst density profile. The spike is of width 3.15 x 10^^ cm at a radius corresponding to jenc = 10, 
with an initial density jump of factor 10^ pext (Rene), followed by a drop of 10~^ pext (Rene)- Left: plot of the forward shock velocity in terms of ^j over radius for 
our numerical simulation and the analytical model of Mesler et al.n2012| . Right: X-ray light curve calculated from our numerical simulation of this setup. 



Using the equations in this section, we are able to set up the simulation grid for optimal results. 

CIRCUMSTELLAR SHELLS 



It has been proposed that circumstellar shells could cause flares in light curves ([ Mesler et al .||2012|). We test this specific 
scenario of the interaction of a blast wave with a narrow shell in which the blast wave becomes non-relativistic. We model this 
scenario as a stellar wind environment followed a large den sity jump of small w idth that is terminated by a large density drop 
followed by a homogeneous medium. We model the setup of [Mesler et al.| ( [2QT2] ) but exclude the wind termination shock for we 
have analyzed the effect of a wind termination shock in this paper and have concluded that will not cause a flare. Although we 
do not model the wind reverse shock in this scenario, the reverse shock acts on the mass ejected by the stellar wind. The Lorentz 
factor at the encounter with the strong density jump, according to the swept up mass approximation, is therefor identical with or 
without accounting for the wind reverse shock. Our stellar wind environment setup is: 



10 ergs, mwind = 10 M^yr , v^ind = 10 km/s, and p^ 



Wef- 



^wind 



47r/?2 



.y^wind 



The blast wave then encounters the spike at ^yenc = 10, and the density of the spike is 10^ pext (Rene), where pext is that of Equation^ 
The spike contains the amount of mas s Mspjke = 0. IM q , and is followed by a sudden lower density region of 10~^ pext(Rene)- This 



is essentially the same setup as that of Mesler et al. !(2012) that resulted in light curve flares in their analysis. 

We numerically simulate this setup with a resolution in time of roughly 30 snapshots during the encounter to ensure the behavior 
of the blast wave wit hin the spike is resol ved, an d no flare is observ ed (Right Figure in Figure [12]). The discrepancy between 
our results and that of [Mesler et al. ( [2012[ ) is that [Mesler et al.[ ( [2012 ) models only swept up mass and does not account for the 
radial structure of the blast wave. With this approximation, the mass swept up after the encounter instantaneously mixes with the 
mass swept up prior to the encounter (which are actually separated by the contact discontinuity) resulting in an overestimation 
of the total mass radiating with the post-encounter Lorentz factor In addition to the over-estimation of the mass in the forward 
shock region, this method also overest imates the Lorentz f actor during the encounter To show this, we analytically modeled 
the Lorentz factor evolution as done in Mesler et al. ( [2012[ ) (who follow Pe'er 2012), and compared it to the numerical results 



for the same setup. This is shown in the Left Figure m Figure 12 This figure depicts the overestimation of the Lorentz factor 



from the assumptions of Mesler et al. (2012), and the underestimation of the Lorentz factor in the homogeneous medium. The 
Lorentz factor jumps up immediately after exiting the spike, but there is still no flare from this because initially, only negligible 
mass is radiating with this larger Lorentz factor The light curve from our numerical simulation is plotted in the Right Figure 
in Figure 12 and shows a steepening in post-encounter slope but does not show a flare. No flare is observed because the fluid 
behind the contact discontinuity is still traveling at the pre-encounter Lorentz factor immediately after the FS hits spike. The 
FS Lorentz factor immediately drops, resulting in it only minimally contributing to the total flux observed. In the swept up 
mass approximation, the FS Lorentz factor gradually lowers resulting in it still contributing noticeable amounts to the light curve 
immediately after the encounter This may be the cause of the flares seen in Mesler et al. (2012), however, when calculating the 
emission from a spatially resolved blast wave, while taking light travel times into account from different radii and angles, we find 
that no observable flare is produced. 

The light curve shown here is calculated at the X-ray frequency of 5 x 10^^ Hz without the contributions from electron cooling 
or absorption. However, we did calculate light curves with the contributions from electron cooling and absorption, and saw no 
qualitative differences. We also considered light curves for different frequencies such as optical and radio. Light curves calculated 
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at optical frequencies exhibited the same behavior as the X-ray Hght curves. Below the synchrotron self- absorption break a steep 
drop was observed following the onset of the shell encounter and a steep rise when the blast wave emerged from the shell. 
However, the (thin) emitting region of the blast wave in the optically thick regime was not always fully resolved numerically, 
leading to noise in the light curves. Regardless, our conclusions regarding X-ray light curves, which are the main focus of our 
study, remain unaffected. 

We believe that any model where the blast wave itself is not resolved is fundamentally unsuited to describe features occurring 
on time scales AR/c. In the non-relativistic case, the relevant light crossing time is of the order of the blast wave radius, R/c, 
due to both the absence of relativistic beaming of the emission and the width of the blast wave AR becoming a sizable fraction 
of the radius. 
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